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Abstract — In this paper, the flow dynamics around two rotating circular cylinders arranged in two 
different configurations, in turbulent and laminar regimes are analyzed through the Immersed Boundary 
Methodology. The Large Eddy Simulation with dynamic Smagorinsky sub-grid scale is also used here. For 
the side by side and in tandem configurations, the simulations are performed, for Re 100, specific rotation 
0.5 and spacing ratio varying between 1.5 and 5.0. For the first arrangement, simulations are also carried 
out, for Re and spacing ratio constant, and the specific rotation varying between 0 and 2. For in tandem 
configuration, the Re ranges from 200 to 10 s and the others parameters are maintained constant. The 
results showed that the cylinder arrangements, in addition to the mentioned parameters, play an important 
role in the mechanism of the vortex shedding, as well as in the wake pattern and in the fluctuations of drag 
and lift coefficients. The rotation mechanism inhibits the vortex shedding process for different values of 
spacing ratio depending on the arrangement of the cylinders. Vorticity contours, time evolutions of fluid 
dynamics coefficients as well as pressure distribution along the cylinder are presented. 

Keywords — Immersed Boundary Methodology, Reynolds number, rotating cylinder, spacing ratio, 
specific rotation. 


I. INTRODUCTION 

The fluid flow with or without heat transfer is 
practically involved in all energy production processes, in 
the environmental phenomena, in the thermal equipment, 
in the aeronautical and aerospace engineering, in the 
reactor engineering, in the bioengineering and so on. In 
projects for automobiles and their components, in projects 
for rotating machines, in the sizing of anti-fire systems in 
large spaces, in the project of refrigeration equipment, in 
the forecast of pollution caused by chimneys in the 
atmosphere and by the discharge of pollutants into rivers, 
lakes and soils, in solving numerous multiphase flow 
problems found in the petroleum industry, etc., numerical 
simulation is of great importance [1], 

Thus, the Computational Fluid Dynamics (CFD) has 
been widely used to study flows around bluff bodies, in 
several analyzes, as mentioned below: Silva et al. [2] 
investigated the flow dynamics around a pair of cylinders 
in different geometric arrangements for Re = 72000; Silva 
et al. [3] analyzed the rotation oscillation effect on the flow 
characteristics at low Reynolds number; Lima et al. [4] 
investigated vortex induced vibration in a circular cylinder 
mounted viscoelastically transverse to flow at moderate 
Reynolds; Dienstmann et al. [5] formulated a simplified 


model for consolidating poromechanical analyses, induced 
by a rigid cylinder in rotating motion immersed in a porous 
medium; Xia et al. [6], investigated numerically through 
the BGK lattice Boltzmann method, the Poiseuille flow of 
the fluid power law on a free rotating cylinder located 
eccentrically in a two-dimensional channel; Wang et al. 
[7], experimentally studied the flow around a circular 
cylinder turning retrograde near the boundary layer of a 
turbulent wall to Re = 10000, using the particle image 
velocimetry technique; Maurya et al. [8] used a rotating 
cylinder arranged in a rectangular T-shaped channel to 
investigate the heat transfer and momentum characteristics 
of a Bingham plastic fluid; Musmar et al. [9] investigated 
experimentally, the effect of rotation increment of a 
cylinder on the wear rate of a cylinder-piston system. The 
engine speed, load and the rotation angle of the cylinder 
were the main parameters; Zou et al. [10] investigated 
numerically the flow induced vibration of a rotating 
circular cylinder for reduced velocity in the intervals 3.0 < 
U* < 14, with the cylinder free to move in the flow 
direction and transverse. Bouchon et al. [11] presented a 
cut-cell method for the simulation of two-dimensional 
flows past obstacles. They showed the results of flows 
around an impulsively started circular cylinder for 
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Reynolds number 1000 and 3000 and also, results of flows 
around a moving rigid cylinder at Re = 800. The semi- 
staggered grid technique was employed by Shu et al. [12] 
to solve several classical testing cases as flow in a channel 


2.1 Mathematical formulation for the fluid 

The filtered Navier-Stokes equations and the continuity 
equation for viscous, incompressible and Newtonian flow, 
can be presented respectively, in tensorial form, as follow: 


lid-driven flow. 

du i ^{u i u j / 

1 dp | 


du , du, 

It is worth mentioning here, the Immersed Boundary 

dt dXj 

Pdx,, 

V cf 

dx. dx. 

\ j 1 J 


+ /, 


Methodology (IBM), which has increasingly gained 
prominence in the study of flow of simple, complex, 
mobile and even deformable geometries, without need for 
rewrapping. Initially proposed by Peskin [13], it has been 
used and improved for different applications, such as in 
hydrodynamic magnet flow (MHD) with complex and 
mobile boundary problems [14]; in simulation of an elastic 
vesicle close to a wall in shear flow using the modified FI 
method [15]; in the mass diffusion and convection 
modeling through porous membranes under large 
deformations [16]; in problems involving heat transfer 
[17]; in the study of linear elliptic value problems on 
arbitrary domains with a Dirichlet boundary condition was 
proposed an augmented Immersed Boundary Method [18]; 
in the investigation of the flow physics behind dynamical 
transitions that occur in the flow field around a plunging 
foil [19]; among others. 


^ = 0 
dx, 


( 1 ) 


( 2 ) 


Where p [Kg/ in 1 and v ef [m 2 / s] are respectively, 

the specific mass and the effective viscosity, properties that 
characterize the fluid. The interest variables are 
represented by: u i \m / ,v] is i-th component of the filtered 

velocity; p |^A//h 2 J is the filtered pressure and 

f |^A/m 3 J is the i-th component of the Eulerian force, 

that acts on the fluid. The source term of Eulerian force 
only exists in the Eulerian points coincident or close of the 
Lagrangian grid, being null for the other points of the 
domain. It is calculated by distribution of the components 

of the Lagrangian interfacial force vector F ( x k ) [ N ], 


In the present work, numerical simulations are 
performed for the analysis of two-dimensional, 
incompressible and isothermal flow around a pair of 
cylinders arranged side by side and in tandem. The 
influence of the spacing ratio is investigated for the two 
cylinders configurations, and the influence of the Reynolds 
number is analyzed for the aligned cylinders. Also, the 
rotating movement effect for the side by side case is 
investigated. The vorticity fields, the time evolution of the 
fluid dynamics coefficients, as well as the distribution of 
the pressure coefficient, are presented. The Immersed 
Boundary Methodology, together with the Virtual Physical 
Model is used (VPM) [20], 


using a distribution function, as proposed by Peskin and 
McQueen [21]. The equation for Lagrangian force can be 
expressed as: 

F(x k ,t) = p SV ^ + P V.(v(x k ,t)V(x k ,t))- 

^ p * 

, r * _ (3) 

V.| p(vv (x k ,t) + < / T V(x k ,t)) | + V P (x k ,t) 

In (3) F a [A] is the acceleration force, F i [a] is the 

inertial force, F v [A] is the viscous force and F p [A] is 
the pressure force. 


II. THE MATHEMATICAL MODEL 

The IBM [13] together with VPM [20] are used in the 
present work, to simulate two-dimensional flows around 
two rotating circular cylinders arranged side by side and in 
tandem configurations. It is based on the Navier-Stokes 
equations plus a force term that acts so that the fluid 
“perceives” the interface, thus making the exchange of 
information between fluid and solid. In this method, the 
Eulerian grid represents the calculation domain and the 
Lagrangian grid, the immersed interface. These grids are 
geometrically independent and coupled through the force 
term. 


2.2 Indicator function 

The indicator function, proposed by Unverdi and 
Trygvason [22] is also used in the present work. It is an 
interface tracking method, where the function is calculated 
in whole domain or in part of it, with the attribution of the 
unit value for internal points to the interface, zero for the 
external points and values between zero and unit for the 
transition points, that is, on the interface. This function is 
based on a function G(x,f) and can be expressed by: 

V/(jf,r) = G(x,r) ( 4 ) 

The term on the right side of (4) is given by: 
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G(x,t) = Y J D ij {x-x k )n(x k )AS(x k ) ( 5 ) 

k 

In (5) h(x k ) is the normal vector of the interface. 

Applying the divergent operator in (4), the Laplacian of 
the indicator function is obtained: 

V 2 /(x,f) = VG(x,f) ( 6 ) 

In this way, after solving (6), the solution of the 
indicator function is obtained, / (x, f) , in all points of the 

domain. For the resolution of the linear system resulting 
from the discretization of (6), the SIP solver is used. 

2.3 Turbulence model 


of p that allow that the values of the velocity components 

u n+ ' and u" 11 , obtained from the respective Navier-Stokes 
equations, satisfy the mass conservation in the instant of 
time n +1. Thus, the estimates of the velocity components, 
u n+1 , obtained with information from the fields u , v , p 
and / of the preceding time, are: 


~n+l 

ll —u 


At 


— + 


d{u { Uj ) 

dxj 


1 dp" 
p fix,. 


d 

V - 

du t duj 

dxj 

ef 

i dxj fix,. 


(7) 


Turbulence is one of the most challenging problems of 
the modern physical and it is among the most complex and 
fascinating phenomena found in the nature. Due to several 
practical implications in different sectors, the number of 
researches related to the understanding and control of 
turbulent flows has become increasingly greater. An 
alternative for the treatment of flows at high Reynolds 
number can be done by scales separation. This separation 


Making the appropriate arrangements, obtain: 

1 9 {p" + '-p") 


At 


P 


dx , 


( 8 ) 


After the application of the divergent operator — on 

dx i 

both sides of (8), follows: 


can be accomplished through the decomposition process of 

1 

du" +l 

M+l 

OUi 

_ 1 9 

f o 'n+1 "N 

op 

the Navier-Stokes equations into a mean and a floating 

At 

fiX; 

fix,. 

p fix,. 

{ fa, J 


part, as proposed by Reynolds [23], or through a filtering 
process proposed by Smagorinsky [24]. Here, the Large 
Eddy Simulation (LES) is used, for which the filtering 
process is necessary. For the turbulent viscosity 
calculation, used in LES, there are different sub-grid 
models. In the present work, the sub-grid model of 
Smagorinsky is used, which is based on the hypothesis that 
the production of the sub-grid turbulent tensions is equal to 
the dissipation. The turbulent viscosity is given as a 
function of the strain rate and the length scale. The 
Smagorinsky constant used here is 0.18 for isotropic and 
homogeneous turbulence. 


(9) 

On what p " +1 = p" + - p" is the correction pressure. 
The importance of the divergent is due to the fact that it is 
one of the terms in the equation that expresses the mass 
conservation. As it is necessary that the velocity field 
satisfy the continuity equation, the second term on the left 
9 m" +1 

side of (9), — ! — , will be null. Therefore, the equation can 

dx i 

be rewritten as: 


1 dii" 


1 d 2 p ' 


( 10 ) 


At dx :. p dxjdxj 


III. NUMERICAL METHOD 

For the velocity-pressure coupling, the Fractional Step 
Method is used, initially proposed by Chorin [25], It is a 
non-iterative method, where from the fields of u , v, p 

and / of the previous iteration, the velocity fields are 
estimated. With these estimated fields, the pressure 
correction is calculated, through the solution of the linear 
system, by the Modified Strongly Implicit Procedure (MSI) 
developed by Scheneider and Zedan [26], The pressure 
behaves like a Lagrange multiplier in minimization 
problems. The importance of the Poisson equation for the 
pressure correction is that it makes the connection between 
the momentum equations and continuity. It provides values 


Rearranging (10), a Poisson equation for the pressure 
correction is obtained, whose source term is the divergent 
of the estimated velocity: 


d 2 p" +1 

dxjdxj 


p du" +l 
At fix,. 


or V 2 p n+l 


-?-V.u" +l 

At 


( 11 ) 


Thus, the estimated velocity fields are obtained through 
(7) and the pressure correction field, through the resolution 
of the linear system, generated by discretization of (11). 
From (8), the velocity field is updated by solving (12): 


At dp " +1 
p fix,. 


( 12 ) 


For the time advance of the momentum equations, second 
order Adams-Bashforth is used and, for the first iteration, 
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the first order Euler method. Thus, the discretization of the 
components x and y of the estimated velocity is given 
respectively by: 

<j = <J + [!• 5 {-K, + Dh ) - 0,5 (-A£ + D n T \)] - 


-At| -P^j+fl 


(13) 


u n+1 = u n .+At\l,5(-A n \ D" ) - 0, 5 (-A"” 1 + D"7 l )] - 

>,i i,j L V y-i y>j! V y‘<J yo / J 


-At 




(14) 


For the Navier-Stokes equations, the second-order 
centered finite differences in space are used. A displaced 
arrangement is adopted for the velocity components 
because it provides more stability in the pressure-velocity 
coupling. For the discretization of the Lagrangian force, the 
calculation of the existing derivatives, is done through the 
Lagrange polynomials using auxiliary points in the x and 
y directions. Thus, to determine the variables in the 
auxiliary points and in the points of the Lagrangian mesh, 
interpolations are made from the Eulerian mesh variables. 


IV. PROBLEM DESCRIPTION 

The simulations are carried out in a non-uniform mesh 
of 600x300 points, calculation domain 50Dx30D, being D 
the cylinder diameter. Downstream from the cylinder 
center until the end of the domain is higher than 26D. 
According to Prasanth et al. [27], the boundary 
downstream of the cylinder, at a distance equal to or 
greater than 25.5D, has no significant effect on the flow 
and in the cylinder response. On the other hand, Chen and 
Zha [28], after intensive numerical experiments, concluded 
that the solution is not influenced when the downstream 
boundary is located 20D from the cylinder center. 

For the two circular cylinders arranged side by side, 9 = 
90°, the upper cylinder rotates in the clockwise direction 
and the lower cylinder in the counterclockwise direction, 
Fig. (la). For in tandem cylinders case, 9 = 180°, the 
upstream cylinder rotates in the clockwise direction 
whereas the downstream cylinder in the counterclockwise 
direction. Fig. (lb). 


16.5D 


LL 


16.5D 


i (a) 

! (b) 


O O 


L *! 


30D 


50D 

Fig.l: Rotation direction of each cylinder 


The movement is imposed around its own axis through 
the imposition of the velocity components at each 
Lagrangian point. For this, the tangential velocity is 
decomposed. 


V. RESULTS AND DISCUSSION 

The study of the flow around one or more stationary 
cylinders obtained by several authors, is of great 
importance to contribute to a better understanding of the 
flow dynamics when cylinders are in motion. In this 
section, are presented the results of the numerical 
simulations, carried out considering two rotating circular 
cylinders side by side and in tandem, in order to investigate 
the spacing ratio, the specific rotation and the Reynolds 
number effect on the flow characteristics. 

5.1 Influence of the spacing ratio on the flow around a 
pair of cylinders side by side 

Here, the vorticity fields are presented as well as the 
time evolution of the fluid dynamics coefficients and the 
respective mean values, for laminar flow regime. 

5.1.1 Vorticity fields 

Figure 2 shows the vorticity fields for Re = 100, a = 
0.5, spacing ratios 1.2, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0 and 5.0 
and, the Fig. (2g) shows the vorticity contour, obtained by 
Yoon et al. [29]. 

For the spacing ratio 1.2, Fig. (2a), the two cylinders 
behave as one, due to the great proximity between them 
and the rotation movement in opposed directions. In this 
way, shear layers move from the upper and lower cylinders 
and merge to form a single wake downstream. With the 
increase of the spacing ratio. Fig. (2b), the flow behaves 
very different from the one observed in Fig. (2a). 
Downstream, positive and negative vorticity can be 
observed at the bottom and at the top sides, respectively 
from both cylinders; however, the vortices are not detached 
close to the cylinder, remaining as a single elongated 
vortex. It is even reasonable to say that for L/D = 1.5, the 
rotation movement of the two cylinders inhibited the 
vortex shedding process, even with a rotation equal to 0.5, 
which is confirmed in Fig. (3a). In the results presented by 
Xiao-hui et al. [30], for Re = 450 and T/D = 1.11, the 
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suppression of vortex shedding was obtained for rotation 
equal to 1.74. In Fig. (2c), the shear layers on both sides of 
each cylinder are visible, however, one cylinder influences 
the other so that the vortex shedding process until is not yet 
completely independent. For this spacing ratio, the flip- 
flopping pattern is observed. It is also verified by Yoon et 
al. [29], for stationary cylinders and g* = 0.7, Fig. (2g). For 
L/D = 2.5 the shear layers from both cylinders cause in 
phase synchronized vortex shedding in the first 
dimensionless times, approximately 50 and after the vortex 
shedding is out of phase, as shown in Fig. (3j). It is also 
noted that the opposite rotational movements contribute to 
the formation of a single wake more or less elliptical 
downstream of the cylinders. It is worth mentioning that, 
far from the two cylinders, the wake is narrower than close 
to them. This behavior is reported in Zdravkovich [31], for 
side by side stationary cylinders and 1.1-1.2 < T/D < 2-2.2. 
The same does not happen for L/D = 3.0, in which just 
downstream, alternating vortices, in phase, are detached 
from each cylinder, in the first dimensionless times, 
approximately < 200. During this interval of time, the 
vortex shedding process is regular and the wake formed 
downstream of two cylinders is clearly independent. After, 
the vortex shedding process, goes through a transition so 
that two independent wakes no more exist giving way to a 
wake with shape more or less elliptical. The vortices are 
also released out of phase and, along of the wake, several 
vortices pairs of lesser intensity are noted. For L/D > 3.5, a 
symmetrical pattern of anti-phase flow is observed. This 
characteristic is also seen for the case of two stationary 
cylinders side by side, for T/D = 4.0 [32]. However, as the 
cylinders have opposite rotation movement, although it is 
verified through the vorticity fields, the fluctuations 
produced by the vortex shedding process are in phase, as 
shown by the time evolution of the drag coefficients for 
L/D > 3.5. Here, for these spacing ratios, two vortices 
wakes are also independent as described by many authors, 
for the stationary case. It is also pointed out, that the 
opposite movement of the two cylinders leads to 
symmetrical wakes in relation to the central line of the 
flow. It can be concluded from Fig. 2, that the spacing ratio 
as well as the rotation movement, have an important effect 
on the vortex shedding pattern, in the separation of the 
shear layers and consequently on the wake pattern. 


(a) LflD=12 

* 1 

(b) L0D=li 

(o) LD=2 jO 

(d) VD=as 

g* -> 

(e) 1A>=20 

(f) Lt>=3 3 

• • • • 1, 

*"*»•*•• | 

(8) 

* • « • • 

(If) LC=5 jO 

r> * . . . * ( 

.. 


(g) 


Fig. 2: Vorticity fields for Re = 100, a = 0.5 and 1.2 <L/D 
< 5.0. (g) extracted from Yoon et al. [29], for a = 0.0 and 
g* = 0.7 

5.1.2. Time evolution of the fluid dynamics 
coefficients 

Figure 3 shows the time evolution of the drag 
coefficient, for the same data as in Fig. 2. The left column 
corresponds to the lower cylinder and the right column to 
the upper cylinder. It is observed that the fluctuations show 
similar behavior for the two cylinders, for all spacing 
ratios. 

For L/D = 1.2, the fluctuations are periodic and no 
sinusoidal. It is worth mentioning that, although the 
amplitudes of these, for the both cylinders have the same 
tendency, the fluctuations are out of phase. For L/D = 1.5, 
the fluctuations amplitude decreased drastically, keeping it 
approximately uniform throughout the simulation. This 
confirms the great influence of the rotation in the vortex 
shedding process, for this spacing ratio, as shown in Fig. 
(2b). For L/D = 2.0, the fluctuations show an irregular 
behavior over the time, corresponding to the flip-flopping 
pattern. For the spacing ratio equal to 2.5, the fluctuations 
are synchronized in phase for the dimensionless time 
approximately up to 50, as already commented. After, 
“kinks” appear in the fluctuations, which are similar for 
both cylinders, just lagged in time. Fig. (3j). This behavior 
is due to the shape more or less elliptical of the wake, far 
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from the cylinders. Fig. (2d), and due to vortex fusion, that 
occurs during the vortex shedding process. For L/D = 3.0 
the fluctuations for the drag coefficient are periodic, 
sinusoidal and of great amplitudes for the dimensionless 
time < 200. After these times, the fluctuations remain 
regular, but with less amplitude. This behavior is due to the 
change in the wake pattern, throughout the simulation, as 
previously commented. For the spacing ratios > 3.5, the 
drag fluctuations present similar behavior for both 
cylinders, with slight reduction in amplitude as the spacing 
is increased. It is interesting to point out that although the 
vorticity fields show an anti-phase flow pattern. Fig. 2, as 
classified for the case of two side by side stationary 
cylinders, here, as the cylinders are in opposite rotation 
movement, the drag fluctuations for the two cylinders are 
synchronized in phase, as already commented. For a better 
understanding of what is being exposed, a zoom has been 
done, considering L/D = 4.0 for the fluctuations of the both 
cylinders. Fig. (3i). 
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Fig. 3: Time evolution of the drag coefficients, for Re 
= 100, a = 0.5 and several spacing ratios. 


Figure 4 shows the time evolution of the lift coefficient, 
for the same data of the Figs. 2 and 3. The left column 
corresponds to the lower cylinder and the right column to 
the upper cylinder. The behavior presented for all spacing 


ratios is also similar for both cylinders. It is interesting to 
point out that the vortex shedding process and the wake 
pattern have different effects on drag and lift coefficients. 
For spacing ratio equal to 2.5, the lift fluctuations are in 
phase, as shown in Fig. (4i), contrary of drag fluctuations 
that are lagged in time. Fig. (3j). For L/D > 3.5, again 
opposite to what is observed for the drag fluctuations, the 
lift fluctuations are anti-phase, as can be seen by zoom 
made for L/D = 4.0, Fig. (4j). It is worth mentioning that 
for all spacing ratios, the lift obtained for the lower 
cylinder (counterclockwise rotation) is negative and the 
obtained for the upper cylinder (clockwise rotation) is 
positive. This behavior is expected once a time that, the 
velocity at the top of the upper cylinder (same direction of 
the flow) is greater than the bottom (opposite direction of 
the freestream velocity). 


As a result, the pressure at the bottom is greater, 
causing a positive lift force. On the other hand, for the 
lower cylinder, the velocity at the top is in the opposite 
direction of the flow, being, therefore, smaller than that the 
bottom (same direction of the freestream velocity). 
Consequently, the higher pressure at the top, causes a 
negative lift coefficient. 
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Fig. 4: Time evolution of the lift coefficients, for Re = 
100, a = 0.5 and several spacing ratios. 
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Figure 5 shows for better understanding, the mean 
values of the fluid dynamics coefficients for the two 
cylinders, placed side by side. With the exception of what 
is observed for L/D = 1.5, the mean values of the drag 
coefficients for both cylinders, are practically the same. On 
the other hand, the mean values of the lift coefficients have 
opposite behaviors, due to the direction of the cylinder 
rotation. Here, the positive lift generated by the upper 
cylinder with clockwise rotation, and, the negative lift 
generated by the lower cylinder with counterclockwise 
rotation, are better visualized. 



2 - 
1 - 


-1 ■ 
-2 - 


upper cylinder 


Cl 


/" Icwer cylinder 


2 l/d 3 


Fig. 5: Mean values of the drag and the lift 

coefficients as a function of the L/D for side by side 
cylinders. 

5.2 Influence of the specific rotation on the flow around 
a pair of cylinders side by side 

In this section the vorticity fields as well as the mean 
values of the drag and lift coefficients are presented, for Re 
= 100, L/D = 1.2 and 0 < a < 2.0. 

5.2.1 Vorticity fields 

The flow visualization by the vorticity contours is 
shown in Fig. 6, after the flow has reached the permanent 
regime. It can be noted the negative vorticity generated by 
the upper cylinder and the positive vorticity generated by 
the lower cylinder, as a single body. One probable reason 
for this behavior is that the cylinders are rotating in 
opposite direction from each other. The interaction 
between two shear layers originated from each cylinder 
occurs only at downstream of the cylinders, as expected. 
There is a ‘2S’ mode of vortex shedding in the classical 
Von Karman Street. Then, the wake formed behind them is 
symmetrical in relation to the center of the domain. As can 
be seen, the results show that the cylinder rotation has an 
important effect on the vortex shedding process. This 
mechanism is suppressed as the specific rotation increases. 
In this way, the flow reaches a steady state, without 
vortices, at the critical specific rotation equal to 1.3. 
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Fig. 6: Vorticity contours for Re = 100, L/D = 1.2 and 0 < 
a <2.0. 


The mean values of the fluid dynamics coefficients of 
the present study are compared with the numerical results 
of Yoon et al. [29]. Fig. (7a) shows the mean values of the 
lift coefficients for the upper and lower cylinders, and Fig. 
(7b) the mean values of the drag coefficients. The lift 
coefficients are equals for both cylinders, however with 
opposite values. The lower cylinder has negative mean 
values while the upper cylinder presents positive mean 
values, for all analyzed specific rotation. 

It can be observed. Fig. (7a), that the mean values 
decrease with an increase of the specific rotation for the 
lower cylinder, presenting a little increase for a = 2.0. On 
the other hand, the values increase as the specific rotation 
is increased for the upper cylinder. In Fig. (7b) it is noted 
that the mean values of the drag coefficients decrease with 
an increase of the specific rotation. This is coherent once a 
time that, the vortex shedding mechanism is suppressed as 
the specific rotation increases, as shown in Fig. 6. The 
obtained results presented good agreement with the 
numerical results of Yoon et al. [29]. 
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5.3 Influence of the spacing ratio on the flow around a 
pair of cylinders placed in tandem 


Here, the vorticity fields are presented, as well as the 
time evolution of the fluid dynamics coefficients and the 
respective mean values, for laminar flow regime Re = 100, 
a = 0.5 and spacing ratios equals to 1.5, 2.0, 2.5, 3.0, 3.5, 
4.0, 4.5 and 5.0. 


5.3.1 V orticity fields 

The flow visualization by vorticity contours is shown in 
Fig. 8 after the flow has reached the permanent regime. 
Here, the upstream one has a clockwise rotation movement 
and the downstream one has a counterclockwise rotation. 
For spacing ratio from 1.5 to 3.5, Figs. (8a) to (8e), the 
both cylinders behave as a single bluff body and the vortex 
shedding process reduces as the spacing ratio increases. It 
is also noted that the vortices behind the cylinders are more 
elongated, which reduces the drag on the upstream cylinder 
as the spacing increases. On the other hand, as the 
downstream cylinder is surrounded by the shear layers of 
the upstream cylinder, consequently in a region of lower 
pressure, the drag on it is significantly less than that of the 
upstream cylinder. It interesting to point out that, for L/D 
equals to 3.0 and 3.5, the vortex shedding process occurs 
normally for the first times of simulation and decreases as 
time advances. For spacing ratios from 4.0 to 5.0, the 
vortex shedding process of the upstream cylinder becomes 
more independent, as well as that of the downstream 
cylinder. It is also verified that the generated vortices are 
closer to each other, are more robust and have a rounded 
shape, differently from what is observed for smaller 
spacing (vortices more elongated). As a result, the drag on 
two cylinders increased. From the presented results, it is 
clearly confirmed that the influence of the spacing ratio as 
well as the rotation movement in the suppression of the 
vortex shedding mechanism, depends on the cylinders 
arrangement. 
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Fig. 8: Vorticity contours for Re = 100, a = 0.5 and 1.5 < 
L/D < 5.0. 


5.3.2 Time evolution of the fluid dynamics 
coefficients 

Figure 9 shows the time evolution of the drag 
coefficients, for cylinders in tandem, Re = 100, a = 0.5 and 
1.5 < L/D < 5.0. As previously mentioned, the vortex 
shedding process is influenced by the spacing ratio and the 
rotation movement of the cylinders. Consequently, it also 
reflects in the amplitude of the fluctuations of the drag 
coefficients for both cylinders, which reduces as the 
spacing ratio increases up to 3.5, even keeping the constant 
rotation. Thus, for L/D = 3.0 and 3.5, the vortex shedding 
process occurs only in the first times of simulation. As time 
progresses, even for these spacing, due to the influence of 
the rotation movement, the vortex shedding process has 
been suppressed, consequently, the drag fluctuations as 
well. For L/D > 4.0, the amplitude of the fluctuations 
increases significantly, compared to those observed for 
small spacing ratios. This is consistent, once a time that the 
generated vortices are stronger and more rounded, Figs. 
(8f) to (8h). It is also noted that the fluctuations amplitude 
remained practically the same with the increase of the 
spacing ratio, for both cylinders. 
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Fig. 9: Time evolution of the drag coefficients, for 
upstream cylinder (left column) and downstream cylinder 
(right column). 


Figure 10 shows the time evolution of the lift 
coefficient for the same spacing ratios as in Fig. 9. For L/D 
= 3.0 and 3.5, fluctuations for the dimensionless time 
greater than 100 are no longer observed. On the other hand, 
for 3.5 < L/D <3.0 the rotation movement had no effect in 
suppressing the vortex shedding process. It is noted that, 
for the upstream cylinder, the fluctuations amplitude 
decreases as the spacing ratio increases up to 3.5 and 
increases with further increase in L/D. This behavior is 
coherent once a time that, for L/D > 3.5, the vortex 
shedding mechanism occurs normally. The same behavior 
is verified for the downstream cylinder. However, the 
amplitude of the fluctuations is greater. 



Fig. 10: Time evolution of the lift coefficients, for 
upstream cylinder (left column) and downstream cylinder 
(right column). 


It is worth mentioning that, the upstream cylinder has a 
positive lift coefficient, once a time that, has a clockwise 
rotation movement, whereas the downstream cylinder has a 


negative lift coefficient, due to the counterclockwise 
direction movement. 

Figure 11 shows the mean values of the drag and lift 
coefficients for the upstream and downstream cylinders, 
aiming at a better understanding of the previous Figs., (9) 
and (10). It can be seen that the drag exerted on both 
cylinders has opposite effects for L/D up to 3.5, this is, for 
the upstream cylinder, the drag reduces as the spacing 
increases and for the downstream cylinder, the drag 
increases. With further increase in the L/D, there is an 
increase in the mean value of the drag for both cylinders, 
being more significant in the downstream one. The 
explanation for this behavior is well visualized in Fig. 8, in 
which, for L/D = 3.5, the vortex shedding process is 
practically suppressed, due to the rotation movement, while 
for L/D > 3.5, the wake vortices is visualized. The mean 
values of the lift coefficient gradually reduce with 
increasing spacing ratio up to 4.0 for the downstream 
cylinder and then gradually increase as well. For the 
upstream cylinder, this increase is observed from a smaller 
spacing ratio (L/D = 3.0), after which the mean value 
grows and remains practically constant after L/D = 4.0. 



Fig. 11: Mean values of the drag and lift coefficients in 
function of L/D 


It is interesting to note, that the rotation movement of 
the cylinders, produces more significant changes in the 
drag than in the lift, as shown in Fig. 11. 


5.4 Influence of the Reynolds number for in tandem 
cylinders 

Figure 12 shows the time histories of the drag 
coefficient for the upstream cylinder, considering the 
laminar and turbulent flow regime, L/D = 2.0, a = 0.5. It is 
clearly noted that the drag on the cylinder decreases 
continuously as the Re increases. Fig. (12a). For Re < 500, 
the fluctuations are periodic, which corresponds to the 
alternating vortex shedding process, with the formation of 
the classical Von Karman Street downstream of the two 
cylinders. 
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Fig. 12: Time evolution of the drag coefficients, for 
upstream cylinder, 200 <Re< 10 s . 


For Re 1000 and 2000, the fluctuations pattern is 
similar. It is noted the presence of “kinks” in the 
fluctuations along the time evolution, indicating the 
generation of more than one vortex per cycle. The wake 
seems disorganized, with vortices pairs scattered along it 
and, with different size and strength. The “kinks” pattern 
changes for Re = 3000, indicating a greater number of 
weak and strong vortices per cycle. It is interesting to note 
that for Re = 4000, there is a reduction in the mean value 
of the drag coefficient throughout the simulations. For the 
dimensionless time approximately between 40 and 173, the 
mean value is 0.9585, while for the rest of the time it is 
0.9248. In addition, there is a reduction in the amplitude of 
the fluctuations for the mentioned time intervals. As 
increases the Re, the increase in the amplitude in the first 
dimensionless time became more considerable, being 
accompanied by a sharp reduction. In addition, the 
presence of envelopes, in a greater or lesser amount, for the 
cases with Re > 5000, is observed. 

It can be seen from obtained results, that although the 
numerical simulations are two dimensional, with the 
increase in the Re and the preponderance of the inertial 
forces over the viscous forces, three-dimensional 


instabilities develop, turbulence levels begin to appear, 
which influences the vortex shedding process, even though 
the periodicity still remains robust [33]. Such instabilities, 
as well as the totally disorganized vortex street certainly 
reflect in the drag and changes the fluctuations pattern as 
shown in Fig. 12. 


Figure 13 shows the time evolution of the drag 
coefficient for the downstream cylinder. The shear layers at 
the top of the upstream cylinder, in a clockwise rotation 
movement, collide with the downstream cylinder, while 




Fig. 13: Time evolution of the drag coefficients, for 
downstream cylinder, 200 <Re < 10 5 . 


those originating from the bottom, surround the 
downstream cylinder, causing smaller drag and greater 
fluctuations amplitude. 

For the two first simulated Re, the fluctuations pattern 
is similar to that observed for the upstream cylinder. In 
addition, it is clearly noted, that the mean drag is positive 
for Re = 200. With the increase of Re up to 3000, the 
fluctuations pattern is no longer sinusoidal and differs from 
that observed for the upstream cylinder. And, it is noted 
that the mean drag has a negative value. For Re = 4000, 
two different fluctuations pattern are verified, which 
implies that the vortex shedding pattern has changed over 
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the simulations. With further increase in the Reynolds 
number, there is a considerable increase in the amplitude of 
the fluctuations, for the first dimensionless time and after, a 
sharp reduction. This behavior is also observed for the 
upstream cylinder. 


Figure 14 shows the time histories of the lift coefficient 
for the upstream cylinder. Note that the lift coefficient is 
positive for all analyzed Reynolds number. This is 
coherent, once a time that, as the upstream cylinder rotates 
in the clockwise direction, the velocity at the top of the 
cylinder is greater than at the bottom. 




100 Tadhn 


Fig.14: Time evolution of the lift coefficient for the 
upstream cylinder, 200 <Re< 10 s . 


As a result, the pressure at the bottom is greater, 
causing a positive lift force. It is also worth to note, that the 
coefficient reduces as the Re increases, with a slight 
increase from Re 6000 to 8000 and then remains 
approximately constant for the others Re, with small 
variations. It is clearly observed that the fluctuations 
pattern varies along of the simulations, which implies that 
the vortex shedding process has been changed. 

The time evolution of the lift coefficient for the 
downstream cylinder is shown in Fig. 15. Contrary to what 
is observed for the upstream cylinder. Fig. 14, the 


downstream cylinder has a negative lift coefficient for all 
analyzed Re. This is due to the fact that, with the 
counterclockwise rotation movement, the velocity at the 
top of the cylinder is smaller than that of the bottom and, 
consequently, the pressure at the top is greater, resulting in 
a negative lift coefficient. It is also observed that the 
fluctuations amplitude for the downstream cylinder are 
greater than for the upstream one, with a behavior 
approximately similar to that observed in Fig. 14. 




Fig. 15: Time evolution of the lift coefficient for the 
downstream cylinder, 200 <Re < 10 5 . 


5.5 Influence of the Reynolds number in the pressure 
coefficient 

The influence of the Reynolds number can also be 
analyzed through the mean values of the pressure 
coefficient around the surface of the rotating cylinder. The 
same Re values of the Fig. 15 are simulated. However, for 
better visualization, only the pressure distribution for five 
Re values is plotted. In Fig. (16a), are shown the pressure 
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distribution along of the upstream cylinder. It can be seen 
that, at the point of minimum velocity (stagnation point), 9 
= 0°, the local pressure is maximum and consequently the 
mean pressure coefficient also, Cp = 1. At the bottom of 
the cylinder (9 ~ 90°), it appears that the values for the 
mean pressure coefficient are higher than those 
corresponding to the top side, 9 ~ 270°. This is coherent, 
once a time that, the flow is accelerated at the top and 
decelerated at the bottom. With the increase of Re up to 
3000, a reduction in the pressure coefficient can be seen at 
the bottom side of the cylinder, and after this Re, an 
increase is verified, with the mean values being between 
those obtained for Re 200 and 3000. It is interesting to 
observe that for Re > 1000 there is an inflection in the 
curve to the angle approximately equal to 90°. On the other 
hand, at the top of the cylinder, the coefficient shows 
variations with the increase in the Re, with the mean values 
between those obtained for Re 500 and 50000. In addition, 
it is noted that the curve moves slightly to the right. 



Fig. 16: Mean pressure coefficient in function of 
the angle 

The mean values of the pressure coefficient for the 
downstream cylinder are shown in Fig. (16b). It can be 
seen that, at 9 = 0°, the local pressure is minimum and 
consequently the mean pressure coefficient also, Cp ~ -1.0. 
This is due to the fact that the shear layers coming from the 
upper side of the upstream cylinder collide with the front 
part of the downstream one, reducing the local pressure 
and consequently the pressure coefficient. Contrary to what 
is observed for the upstream cylinder, it is noted that the 
lower side of the cylinder, 9 ~ 90°, has smaller values than 
the upper side, 9 ~ 270°, for the pressure coefficient. This 
behavior is coherent, once a time that, the downstream 
cylinder has a counterclockwise rotation movement. Thus, 
the flow is accelerated on the lower side (smaller local 
pressure) and decelerated on the upper side (higher local 
pressure). It is important to note that the stagnation point, 
this is, point of minimum velocity and maximum pressure, 
has changed its position in relation to observed for the 
upstream one. It is verified a displacement of it towards the 
top of the cylinder, as shown in Fig. 17, by visualization of 
the pressure fields, for Reynolds number 2000 and 50000. 


Here, for Re = 200, it is noted that the maximum mean 
value of the pressure coefficient, approximately 0.68, is 
obtained for the angle equal to 310°, less than the 
maximum value obtained for the upstream cylinder, Cp = 
1.0. It is also noted that the Cp presents increasing values, 
for the angles of approximately 20° up to 310° and then it 
reduces in an accelerated way until 360°. Along the 
cylinder, it is verified that there are points at which growth 
is slower than at others. With the increase of Re to 500, the 
flow has a similar behavior, however for the angles 20° up 
to 310°, the mean values showed a reduction and in the 
range of 300° up to 340° an increase, remaining below 1.0, 
approximately 0.85. The maximum value of Cp (0.967) is 
verified for Re = 3000, with the curve moving to the right 
(from 310° up to 316°, approximately). For the others 
simulated Re, the mean values are between those obtained 
for Re = 200 and 3000. 
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Fig. 17: Pressure fields for L/D = 2.0, a= 0.5. 


VI. CONCLUDING REMARKS 

Two-dimensional numerical simulations are carried out 
in order to investigate the flow dynamics around two 
rotating circular cylinders arranged side by side and in 
tandem. The results showed that, beside the influence of 
parameters such as Reynolds number, spacing ratio and 
specific rotation, on the flow behavior, the arrangement of 
the cylinders along the flow is also important. This is 
verified, through the analysis of the vortex shedding 
mechanism for both cases, considering a = 0.5 and Re = 
100. For the side by side case, the suppression of the vortex 
shedding process is obtained for L/D = 1.5. On the other 
hand, for the in tandem case, the vortex shedding gradually 
decreased with the increase in the spacing, and the vortex 
suppression is obtained for L/D = 3.0 and 3.5 after a 
dimensionless time of approximately 100. Thus, it is 
showed the great importance of the cylinder arrangements 
combined with the spacing ratio on the vortex shedding 
process, separation of the shear layers and consequently on 
the wake pattern. It is also verified that, keeping the 
spacing ratio and Re constant, the flow dynamics is 
strongly influenced by the specific rotation. In the present 
work, the critical specific rotation, for which the vortex 
shedding mechanism is suppressed is 1.3. 

It is also found, for the aligned cylinder case, that the 
Re number has a strong influence on the flow 
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characteristics, so that the wake seems disorganized, with 
three-dimensional instabilities and vortices with different 
sizes and forces. In addition, these characteristics reflect in 
the drag and lift force, as well as change the fluctuations 
pattern. Another great influence of cylinder arrangement is 
relating to the stagnation point of the downstream cylinder 
that displace towards the top of the cylinder, whereas in the 
upstream cylinder this point is located at front part, 0= 0°. 

Further studies will be performed in order to analyze 
the flow dynamics around two or more circular cylinders 
with staggered arrangements and also submitted to others 
types of movement such as rotational oscillation and one 
degree of freedom. 
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